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INTRODUCTION 


Recently,  there  has  been  a  revival  of  interest  in  the  dynamics  of  a 
single  bubble  driven  into  pulsation  by  an  applied  acoustic  field.  A 
technique  was  needed  which  would  allow  the  close  examination  of  an  indi¬ 
vidual,  undistorted  bubble.  Methods  which  allow  the  trapping  of  a  bubble 
in  a  sound  field  by  some  mechanical  constraint  do  exist,  but  these  tech¬ 
niques  distort  the  bubble  and  introduce  nearby  boundaries. 

The  purpose  of  this  project  was  to  develop  a  technique  to  study  the 
nonlinear  oscillations  of  bubbles  free  from  the  effects  of  mechanical 
constraints.  The  technique  to  be  described  here  has  been  called  acoustic 
levitation.  A  splendid  isolation  for  the  individual  bubble,  removing  it 
by  over  500  bubble  diameters  from  the  nearest  solid  surface  can  be  obtained 
by  this  technique.  It  also  allows  the  individual  bubble  to  be  viewed  with 
little  difficulty  using  a  small  magnification  microscope,  which  in  turn 
enables  one  to  make  a  precise  determination  of  its  position  and  any 
changes  in  that  position.  Using  this  technique,  the  bubble  may  be 
accurately  positioned  to  within  one-half  of  its  diameter  in  the  sound 
field. 

It  was  desired  that  the  bubble's  size  should  vary  over  an  order  of 
magnitude  for  the  range  of  interest  in  this  study,  and  consequently  a 
method  was  needed  for  measuring  this  range  of  sizes.  Using  the  acous- 


tic  levitation  technique,  the  bubble's  size  can  be  determined  by  measur¬ 
ing  its  terminal  velocity  as  it  rises  freely  through  the  liquid  and 
applying  an  appropriate  drag  law.  Accordingly,  a  wide  range  of  bubble 
"rise  distances"  was  needed.  With  this  technique,  the  range  of  rise 
distances  over  which  the  bubble's  size  could  be  determined  varied  from 
less  than  1  mm  to  nearly  16  mm,  enabling  measurements  of  the  bubble 
radius  to  be  made  from  about  10  pm  to  over  200  pm. 

Acoustic  levitation  occurs  when  an  acoustic  "radiation-pressure" 
force  is  used  to  balance  the  gravitational  or  buoyancy  force  whi^h  tends 
to  force  a  bubble  toward  the  surface.  The  bubble  can  be  trapped  over  a 
range  of  positions  near  a  pressure  antinode  as  a  result  of  the  balancing 
of  these  two  forces. 

The  acoustic  levitation  technique  was  essentially  invented  by  Eller,1 
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and  was  developed  in  great  detail  for  use  with  liquids  by  Apfel.  Crum  ’ 
has  refined  the  technique  and  used  it  extensively  in  the  study  of  gas 
bubbles  in  liquids.  As  an  interesting  sidelight,  the  Jet  Propulsion 
Laboratory^  plans  to  use  acoustic  levitation  of  liquid  metals  in  air 
in  a  future  space  shuttle  mission. 

By  levitating  an  air  bubble  near  the  pressure  antinode  of  the 
acoustic  sound  wave,  it  was  possible  to  observe  small  changes  in  its 
position  as  a  function  of  radius  as  the  bubble  grew  or  dissolved.  These 
measurements  were  then  used  to  investigate  the  nonlinear  oscillations  of 
the  bubble  over  a  range  of  sizes.  The  bubbles  were  studied  in  two  liq¬ 
uids:  isopropyl  alcohol,  and  a  mixture  of  distilled  water  and  glycerine 

(33-1/3%  glycerine  by  volume).  The  temperature  of  the  liquid  was  also 


EXPERIMENTAL  ARRANGEMENT 


The  nonlinear  oscillations  of  bubbles  in  isopropyl  alcohol  and  in 
a  mixture  of  glycerine  and  water  were  studied  by  use  of  an  acoustic 
levitation  technique.  A  detailed  description  of  the  technique  will  be 
given  in  the  experimental  methods  section. 

A  Hewlett  Packard  200  CD,  wide-range  audio  oscillator  was  used  to 
generate  a  signal  in  the  kilohertz  range.  A  Tektronix  545B  oscilloscope 
was  used  to  tune  the  oscillator  to  that  frequency  that  gave  a  maximum 
in  the  acoustic  pressure  amplitude  measured  by  a  pill  transducer  on  the 
side  of  the  levitation  cell.  This  relative  maximum  indicated  that  one 
of  the  stationary  wave  resonances  of  the  cell  was  being  excited.  The 
signal  to  the  cell  was  amplified  by  a  Marantz  8B  power  amplifier  before 
being  introduced  into  two  matched  hollow  cylindrical  transducers  which 
were  part  of  the  levitation  cell.  This  device  consisted  of  the  two 
cylindrical  transducers  separated  by  a  hollow  glass  cylinder.  The  trans¬ 
ducers  were  50  mm  in  height  and  had  a  radius  of  7.7  cm.  The  glass 
cylinder  was  48  mm  in  height  and  had  walls  2.5  mm  thick.  Its  outer 
diameter  was  6.9  cm.  These  three  parts  were  cemented  together  with  a 
thin  glass  bottom  which  acted  as  a  flexible  pressure  release  diaphram. 

The  top  was  covered  but  not  sealed.  A  photograph  of  the  cell  is  shown 
in  Fig.  1. 

The  cell  held  approximately  400  ml  of  liquid.  A  small  pill 
transducer  was  attached  to  the  outer  glass  wall  of  the  cell  to  measure 
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the  acoustic  pressure  amplitude.  The  signal  was  amplified  by  a  Hewlett 
Packard  465A  amplifier  and  its  value  measured  with  a  Hewlett  Packard 
3466A  multimeter.  A  Tektronix  545B  oscilloscope  was  used  to  view  the 
signal  of  the  acoustic  pressure  generated  within  the  cell,  to  insure 
that  a  pure  acoustic  mode  was  generated.  A  block  diagram  and  sketch 
of  the  experimental  arrangement  is  shown  in  Fig.  2. 

A  small,  calibrated  probe  hydrophone  was  used  to  measure  the 
spatial  variation  of  the  stationary  acoustic  wave  along  the  axis  of 
the  cylinder  (the  acoustic  pressure  gradient  is  an  important  parameter 
in  evaluating  the  data).  This  was  done  prior  to  data  collection.  The 
pill  transducer  was  used  in  place  of  the  probe  hydrophone  during  data 
collection.  When  the  bubble  radius  was  to  be  measured,  a  toggle  switch 
was  used  to  turn  off  the  audio  signal  and  start  an  Ortec  873  counter. 

The  jwitch  was  held  in  the  "off"  position  until  the  bubble  was  allowed 
to  rise  a  certain  distance  while  viewed  through  the  cathetometer .  The 
bubble  was  made  visible  by  the  use  of  a  Bausch  and  Lomb  light  source 
which  was  positioned  above  the  cell  and  directed  so  that  light  would 
scatter  off  the  bubble.  When  the  switch  was  returned  to  the  "on" 
position,  the  bubble  returned  to  its  stable  position  near  the  pressure 
antinode. 

The  cell  was  isolated  from  the  room  environment  by  use  of  an  outer 
container,  a  temperature  control  cylinder,  which  fit  around  the  cell 
and  its  leveling  support.  A  Polyscience  Poly temp  90  temperature  control 
system  was  used  to  pump  antifreeze  through  this  cylinder,  maintaining  a 
constant  temperature  in  the  cell.  The  temperature  of  the  liquid  in  the 
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Block  diagram  and  sketch  of  the  experimental  arrangement. 


cell  was  measured  by  the  insertion  of  a  standard  thermometer  into  the 
liquid  through  a  hole  in  the  cell's  cover.  Photographs  of  this  system 
are  shown  in  Fig.  3  and  Fig.  4.  Fig.  3  shows  the  cell  isolation  chamber, 
the  cathetometer  used  to  view  the  bubble,  and  the  illuminator  used  to 
enhance  the  visibility  of  the  bubble.  Fig.  4  shows  some  of  the  electronic 
apparatus  used. 

Standard  methods  to  adjust  the  gas  content  of  the  liquid  were  used. 

An  aspirator  or  a  vacuum  pump  was  used  to  reduce  the  dissolved  gas  con¬ 
tent  of  the  liquid  when  needed.  The  cell  was  set  on  a  marble  slab  and  on 
a  stand  with  adjustable  legs  that  enabled  it  to  be  leveled.  The  stand 
at  the  base  of  the  cell  was  surrounded  by  lead  shot,  between  the 
temperature  control  container  and  the  cell's  base.  The  level  of  liquid 
was  monitored  periodically  in  order  to  maintain  a  sound  field  that  did 
not  change  with  time. 

Because  dirt  particles  present  in  water  could  (and  often  would) 
collect  on  the  surface  of  the  gas  bubble,  it  was  necessary  to  filter 
the  water  before  using.  The  glycerine  and  isopropyl  alcohol  were  clean 
enough  to  be  used  without  treatment.  Batches  of  the  glycerine/water 
mixture  were  made  up  in  large  quantities  and  replacements  taken  from 
this  large  supply  so  that  liquid-property  consistency  was  maintained 


over  several  runs. 


co 

co 

co 


To  collect  a  typical  set  of  data  for  the  analysis  of  the  nonlinear 
oscillations  of  a  bubble,  there  were  several  steps  to  follow.  A  mixture 


of  liquid  was  first  prepared  and  if  necessary  filtered  to  remove  impurities. 
The  liquid's  level  was  selected  so  that  the  pressure  antinode  for  a 
frequency  near  22  kHz  was  underneath  the  transparent  glass  cylinder.  If 
the  liquid  level  was  too  high  or  too  low,  the  bubble  would  assume  an 
equilibrium  position  that  permitted  only  a  limited  length  of  view  as  it 
rose  during  the  measurement  of  its  radius.  The  level  of  the  liquid  was 
adjusted  and  the  correct  frequency  for  that  level  was  determined  for  each 
of  the  different  liquids  used.  However,  once  the  level  was  determined, 
a  small  measuring  device  was  constructed  that  permitted  the  level  to  be 
repeated  after  liquid  replacement. 

Some  data  was  taken  as  a  function  of  temperature;  in  these  cases, 
the  liquid  was  brought  to  the  desired  temperature  before  it  was  intro¬ 
duced  into  the  cell  and  the  temperature  control  system  was  used  to 
maintain  that  temperature.  When  the  temperatures  were  relatively  high 
(50°  C)  or  low  (10°  C) ,  there  was  great  difficulty  in  collecting  data. 
Convection  currents  would  often  be  present  and  would  affect  the  balance 
of  forces  on  the  bubble,  moving  it  in  a  horizontal  or  vertical  direction. 
Even  when  data  was  taken  near  room  temperature  (20°  C) ,  the  temperature 
control  system  was  used  to  help  avoid  any  changes  in  temperature  of  the 
liquid  over  a  data  collection  period.  The  temperature  also  affected 
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the  amount  of  dissolved  gas  in  the  liquid.  If  the  liquid  had  an  excess 
of  gas,  bubbles  would  grow  at  a  much  more  rapid  rate  than  was  desired. 

At  other  times,  the  liquid  would  have  such  a  low  level  of  gas  that 
bubbles  were  hard  to  introduce  and  when  one  was  levitated,  it  would 
quickly  dissolve.  The  requirement  that  the  dissolved  gas  concentration 
be  at  a  particular  value  would  necessitate  the  use  of  one  of  various 
systems  constructed  to  change  the  liquid's  gas  content,  which  in  turn 
would  affect  the  temperature,  which  would  affect  the  gas  content,  etc. 

In  general,  obtaining  the  optimal  liquid  conditions  for  data  acquisition 
was  one  of  many  hours  of  painstaking  trial  and  error  adjustments. 

Another  effect  of  the  temperature  variation  was  that,  at  high  temperatures, 
there  was  some  evaporation  which  caused  the  surface  to  change  after  a 
few  hours,  so  data  had  to  be  collected  and  levels  adjusted  as  quickly 
as  possible. 

Once  the  gas  content  was  close  to  a  level  which vould  allow  the 
slow  growth  or  dissolution  of  the  bubble,  the  frequency  was  constantly 
monitored  to  insure  that  it  was  at  resonance  for  the  particular  standing 
wave  mode  throughout  the  data  collection  period. 

Bubbles  were  introduced  in  various  manners.  One  method  consisted 


of  increasing  the  pressure  until  stable  cavitation  produced  a  small 
bubble.  On  occasion,  a  hypodermic  needle  was  used.  The  introduction 
of  objects  into  the  liquid  required  care  to  avoid  disturbing  the  fluid 
and  changing  its  level.  Once  a  bubble  was  levitated,  care  was  taken 
to  insure  that  only  one  bubble  was  present.  At  times,  several  bubbles 
could  be  levitated  near  the  antinode  and  they  would  then  combine  into  a 
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larger  bubble  due  to  the  secondary  or  mutual  Bjerknes  Force.  Also, 
bubbles  would  appear  on  the  sides  of  the  cell  or  at  the  bottom  or  near 
other  antinodes. 

With  the  individual  bubble  levitated  in  the  stable  liquid,  it 
was  then  necessary  to  collect  data  as  the  size  of  the  bubble  changed. 

It  was  found  that  a  large  bubble,  near  the  resonance  radius,  would 
slowly  dissolve  if  the  gas  content  was  slightly  less  than  saturation. 

A  bubble  which  was  small  would  grow  slowly  to  resonance  size  through 
rectified  diffusion  if  the  liquid  had  the  correct  gas  concentration. 

When  the  gas  conditions  were  properly  set,  it  became  possible  to  obtain 
a  bubble  of  any  desired  size,  and  data  could  be  collected  over  a  complete 
range  of  sizes.  The  initial  size  of  the  bubble  selected  was  determined 
by  visual  observation,  requiring  some  experience,  but  a  nomogram  showing 
rise-times  over  a  unit  distance  for  bubbles  of  various  sizes  was  also 
used. 

With  the  frequency  of  the  standing  wave  and  the  temperature  of  the 
liquid  known  and  with  a  bubble  of  acceptable  size  levitated  in  the 
liquid,  a  data  set  was  taken.  The  location  of  the  bubble  was  recorded 
by  using  the  fiducial  lines  in  the  cathetometer .  The  bubble  was  centered 
on  one  of  the  fiducial  lines,  usually  1,  and  a  reading  was  taken  from 
the  cathetometer  scale,  giving  the  position  of  the  bubble  in  the  sound 
field.  The  sound  field  was  then  turned  off  with  a  toggle  switch,  allowing 
the  bubble  to  move  toward  the  surface.  The  switch  also  started  a  timer 
which  measured  the  time  that  the  sound  was  off.  As  viewed  through  the 
microscope,  the  bubble  was  allowed  to  move  past  a  number  of  fiducial 
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lines  before  the  sound  was  turned  on  again  to  return  the  bubble  back  to 
a  position  near  the  antinode.  The  distance  the  bubble  was  allowed  to 
move  was  determined  by  various  considerations  but  usually  the  size  of 
the  bubble  was  the  most  important.  A  large  bubble  would  move 
more  rapidly  to  the  surface  so  a  long  distance  was  required  in  order 
to  reduce  errors  in  the  rise-time  measurement;  a  rise-time  of  1.5 
seconds  was  usually  the  minimum  value  selected.  The  position  of  the 
bubble  and  the  time  it  took  to  move  a  number  of  fiducial  lines  were 
recorded,  as  well  as  the  voltage,  i.e.,  the  acoustic  pressure  amplitude 
as  measured  by  the  pill  transducer.  In  most  data  sets,  the  acoustic 
pressure  was  held  constant  and  adjustments  were  made  during  data 
collection  only  to  maintain  this  voltage. 

These  steps  were  repeated  for  the  bubble  as  its  size  changed.  An 
attempt  was  made  to  obtain  measurements  of  the  bubble  radius  over  the 
range  from  10%  or  80%  of  the  resonance  radius.  These  measurements  were 
made  on  a  single  bubble  whenever  possible.  A  set  of  approximately  10 
data  points  was  also  made  at  various  voltages  for  small  bubble  sizes 
in  order  to  obtain  a  calibration  curve  for  each  individual  data  set. 

Before  the  data  could  be  analysed  with  the  computer  program,  certain 

physical,  parameters  of  the  liquid  and  gas  system  were  needed.  Most  of 

the  information  could  be  found  in  references  (7)  through  (11),  but  for 

various  temperatures  of  the  glycerine/water  mixture,  the  viscosity  had 

12 

to  be  measured.  This  was  done  using  as  Ostwald  Viscometer  to  measure 
the  time  it  took  the  mixture  to  flow  through  a  given  distance  in  the 
tube  and  comparing  this  to  the  flow  time  of  a  liquid  whose  viscosity 
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was  known.  These  viscosity  measurements  were  made  for  the  range  of 
temperatures  considered.  Other  constants  used  in  the  data  analysis 
included  the  liquid's  sound  speed,  its  surface  tension,  and  its  density. 

A  graph  of  these  values  over  a  range  of  temperatures  is  shown  in  Fig.  5 
for  isopropyl  alcohol,  and  Fig.  6  for  the  glycerine/water  mixture. 

Also  required  in  the  data  analysis  was  the  ratio  of  specific  heats 
of  the  gas,  its  molecular  weight,  and  its  thermal  conductivity.  The 
resonance  radius  for  bubbles  in  the  particular  liquid  was  calculated  by 
another  computer  program  as  a  function  of  some  of  the  various  physical 
parameters.  Various  constants  were  used  that  related  to  the  levitation 
cell,  such  as  the  wavelength  of  the  sound,  and  conversion  factors  to 
relate  the  position  of  the  bubble  in  the  sound  field  to  the  cathetometer 
reading  and  the  fiducial  lines  to  standard  units  of  length.  These  values 
were  obtained  during  a  calibration  run. 

When  all  this  data  (the  17  constants  and  the  measurements  from 
a  data  set)  were  gathered,  they  were  entered  into  a  computer  program  for 
analysis.  From  this  program,  listed  in  Appendix  A  as  Program  1,  the 
following  parameters  of  interest  were  obtained:  a  levitation  number 
(essentially,  the  ratio  of  the  hydrostatic  pressure  gradient  to  the 
acoustic  pressure  gradient),  the  radius  of  the  bubble  for  each  data 
point  entered,  the  polytropic  exponent  of  the  gas  contained  within 
the  bubble,  the  acoustic  pressure  amplitude,  and  the  bubble  damping. 

These  parameters  can  then  be  compared  with  theory  to  test  our  knowledge 


of  bubble  dynamics. 


glycerine /water  mixture. 


THEORY 


In  this  section,  the  theory  relating  to  the  study  of  the  nonlinear 

oscillations  of  bubbles  by  acoustic  levitation  will  be  explained 

briefly.  Since  this  is  not  a  major  portion  of  the  study,  the  reader  is 
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referred  to  previously  published  articles  ’  and  the  computer  programs 
in  Appendix  A-Poly  5  (Program  1) ,  and  Appendix  B-CNSTNT  (Program  2) , 
which  were  used  to  analyze  the  data.  Poly  5  was  used  to  obtain  values 
for  levitation  number  as  well  as  other  information  already  mentioned, 
while  CNSTNT  was  used  to  compute  values  of  the  theoretical  levitation 
number  based  on  a  nonlinear  theory.  Comparisons  of  theory  and  experi¬ 
ment  will  be  made  in  the  data  analysis  section. 

The  primary  Bjerknes  Force,  which  is  the  force  used  to  levitate  the 
bubble  near  the  pressure  antinode,  is  given  by 
FA(r,t)  =  -  <V (t) VP (r , t) > , 

where  the  brackets  on  the  right  hand  side  of  Eq.  (1)  mean  time  average, 
V(t)  is  the  instantaneous  volume  of  the  bubble  and  P(r,t)  is  the  time- 
and  space-varying  pressure  field  which  can  be  approximated  along  the 
axis  of  the  cell  by 

P(r,t)  =  P  -  P.(z)cosRt.  (2) 

CO  A. 

Here  Pm  is  the  ambient  pressure,  Q  is  the  angular  frequency  of  the  sound 
field  and  P  (z)  is  the  space-dependent  amplitude  of  the  stationary  wave, 

A 

with  z  being  the  coordinate  measured  vertically  along  the  axis  of  the 
levitation  cell. 
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If  the  levitation  number  is  large,  it  implies  that  only  a  small 
acoustic  pressure  gradient  is  needed  to  levitate  a  bubble.  This  implies 
that  the  acoustic  force  is  very  efficient  in  stabilizing  the  bubble 
against  gravity.  A  small  levitation  number  implies  the  opposite.  In 
this  study,  the  levitation  number  is  measured  experimentally  and  the 
left-hand  side  of  Eq.  (3)  is  computed  theoretically.  This  indicates 
that  the  results  to  be  shown  in  the  next  section  can  be  considered  as 
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comparisons  between  theoretical  and  experimental  values  of  the  levita¬ 
tion  number. 

The  time-varying  radius  of  the  bubble  is  expressed  by  the  Rayleigh- 

Plesset  equation 

••  3  .2 

p (RR+2  R  )  =  P  -  (P^-PACOsfit)  -  2a/R  -  4pR/R,  (4) 

where  the  surface  tension  is  given  by  a,  the  liquid  viscosity  by  u,  and 
the  instantaneous  pressure  acting  on  the  gas  side  of  the  bubble  inter¬ 
face  is  given  by  P_^.  P^  must  be  determined  by  solving  the  conservation 
equation  relating  the  gas  and  the  liquid  phase.  In  the  case  of 
linear  oscillations,  the  result  of  these  calculations  are  summarized  by 
stating  that  the  thermal  processes  result  in  a  pressure-volume  relation¬ 
ship  of  the  poly tropic  type, 

P±  =  Po(Ro/R)3<, 

where  P  +  2a/R0  is  the  internal  pressure  of  the  bubble  at  equilibrium 
and  k  is  the  polytropic  exponent.  It  has  been  shown'*’^  ^  that  the 
thermal  processes  contribute  a  dissipation  mechanism  which  can  be  repre¬ 
sented  by  modifying  the  value  of  the  liquid  viscosity.  The  acoustic 
radiation  from  the  bubble  also  contributes  to  dissipation,  and  this 
effect  can  also  be  accounted  for  by  modifying  the  value  of  the  viscos¬ 
ity.^  A  description  of  these  thermal  processes  using  a  polytropic 
exponent  has  been  experimentally  verified^  for  small  amplitude  oscilla¬ 


tions,  but  an  extension  of  these  results  for  the  case  of  large-amplitude, 
nonlinear  oscillations  has  never  been  attempted.  Since  the  nonlinear 


oscillations  are  of  concern  in  this  study,  it  is  necessary  to  use  the 
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linear  results  for  k,  and  to  use  the  "effective"  liquid  viscosity, 
although  Eq.  (4)  will  be  treated  nonlinearly, 

For  radial  pulsations  of  moderate  amplitude,  an  approximate  solu¬ 
tion  of  Eq.  (4)  can  be  obtained  by  using 
R  =  R0  [1  +  x  ( t )  ]  , 

18 

and  performing  a  power  series  expansion  of  x.  For  the  linear  case, 
the  lowest-order  approximation  to  the  exact  solution  can  be  obtained 
by  retaining  only  terms  which  are  linear  in  x  and  their  derivatives. 

This  yields 

x  =  ?cos(a)T+6)/[(a)2-u)o2)2+4bV]'5  , 

using  the  following  terminology 

T  =  (Po/p)1*  t/Ro,  W  =  2o/RflPo, 

b  =  2y/R0  (pPo)^,  n  =  VP»’ 

L 

£  =  (l-w)n,  0)  =  R0ft(p/Po)  » 

wo2  =  3k-w,  6  =  tan  [2bu)/(u2-uo2)  ]  • 

The  nondimens ional  resonance  frequency  of  the  bubble  is  given  by 
ido  and  when  it  is  converted  to  dimensional  form,  it  yields  the  expression 
n0  =  (l/RopSpKP^  +  (3<-l)  (o/Ro)  ]  l*.  (5) 

The  effect  of  the  surface  tension's  contribution  is  negligible 
and  k  depends  only  weakly  on  Rg  for  the  range  of  radius  values 
considered;  this  implies  that  R o f2 o  is  approximately  a  constant. 

If  a  resonance  radius  (R^gg)  is  defined  by 

n  =  (I/R^P^HSkP^  (3K-l)(o/RRes)]1'S  (6) 

then  n/fios  ^/^Res"  This  implies  that  for  fixed  sound  frequency  (P.)  ,  the 
ratio  Q/Qo  can  be  varied  by  changing  the  radius  of  the  bubble.  In  the 
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figures  used,  the  levitation  number  is  compared  to  ^/^es  rather  than  to 
fi/flo  since  experimental  measurements  were  made  in  terms  of  variable 
radius  rather  than  variable  frequency. 

3 

To  first  order  in  x,  the  time  averages  in  Eq.  (3)  are  <(R/Ro)  cosf2t> 

s3<xcosfit>,  and  <(R/Ro)3>  =  1*  Thus, 

L  =  4  cosS,  (7) 

e 

where 

B  =  [  (uj2  — to o 2 ) 2  +  4b Vi"15. 


When  the  damping  is  small  and  the  bubble  radii  are  considerably 


lower  than  the  resonance  radius,  Eq.  (7)  becomes 

L  a  (^)  (PA). 

e  2  3k-w  — 


(8) 


CO 

In  this  study  for  the  range  of  bubble  radii  investigated,  w  is  very 
—  2 

small  (w^es  s  10  )  and  k,  the  poly tropic  exponent,  only  weakly  depends 
on  Ro .  This  indicates  -that  the  levitation  number  is  only  weakly  depen¬ 
dent  on  the  bubble  radius  when  R0  is  much  smaller  than  R^eg.  It  has 
13 

been  shown  that  if  the  bubble  is  driven  at  an  amplitude  and  frequency 
such  that  nonlinear  effects  are  significant,  the  levitation  number  has  an 
important  and  measurable  dependence  on  Ro  • 

Harmonic  and  subharmonic  resonances  typical  of  nonlinear  oscilla¬ 
tions  appear  if  higher-order  terms  are  used  in  the  solution  of  Eq.  (4). 
The  solution  x(t)  in  the  harmonic  and  subharmonic  regions  allows  the 
explicit  calculation  of  the  time  averages  in  Eq.  (3).  By  retaining  only 
terms  capable  of  giving  nonzero  contributions  in  the  various  harmonic 
resonance  regions  and  keeping  a  second-order  accuracy,  the  denominator 
of  the  left-hand  side  of  Eq.  (3)  can  be  taken  to  be  simply  1.0. 
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For  the  n  =  2  harmonic  resonance  region,  when  to  =  coo/2,  the  levita- 

3 

tion  number  contains  the  linear  part  2  £Bcos6  plus  other  terms  that 

are  dependent  upon  the  amplitude  and  phase  of  the  resonant  component.  At 

the  n  =  3  harmonic  resonance  region,  where  co  s  u)o/3  it  is  found  that 

Eq.  (3)  reduces  to  the  same  form  found  for  the  linear  approximation 

-f  £Bcos<5.  For  the  main  resonance  n  =  1,  oj  =  u)0,  an  expression  for  the 

L  is  found  in  the  form  of 
e 

Le  =  f  Al  cos<t>!,  (9) 

where  A\  is  the  amplitude  and  is  the  phase. 

Thus,  in  general,  the  levitaion  number  can  be  expressed  as  a  linear 
term  plus  contributions  that  apply  within  individual  resonance  regions 
that  add  additional  terms  to  the  linear  portion  of  the  theory.  These 
harmonic  resonance  contributions  are  extremely  complex  and  won't  be 
given  explicitly  here,  although  the  theoretical  curves  shown  on  the 
graphs  contain  an  evaluation  of  both  the  linear  and  the  nonlinear  effects. 

To  evaluate  some  of  the  equations  used  in  the  study  of  the  nonlinear 
oscillations  which  occur  as  the  bubble  changes  size,  it  is  necessary  to 
have  expressions  for  the  polytropic  exponent  and  the  damping  constant. 

14 

The  results  of  the  linear  theory  for  these  quantities  as  shown  by  Eller 
were  used.  It  has  been  shown^  that  under  normal  conditions,  measurements 
of  the  polytropic  exponent  for  a  variety  of  gases  and  ranges  of  bubble 
sizes  agree  quite  well  with  the  analytical  predictions. 

The  total  damping  constant  b  is  the  sum  of  the  three  principal 
contributions,  i.e.,  viscous,  radiation,  and  thermal.  Thus, 
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b  =  b  +  b  +  b 
v  r  t 


(10) 


where 


bv  =  2u/(pP0K  Ro, 


br  =  15(p/Po)J5(02Ro2/c), 


.h 


b  =  (fi0  Ro/2fi)(p/P0)  dt. 


and 


_  3 (y-1)  [X(slnh  X  4-  sin  X)  -  2(cosh  X  -  cos  X)  ] , 
t  X  (cosh  X  -  cos  X)  +  3  (y-l)X(sinh  X  -  sin  X) 


where  c  is  the  sound  speed  in  the  liquid,  and  y  is  the  ratio  of  specific 


heats  of  the  gas  contained  within  the  bubble.  Also  X  =  R0(2ft/Di)2  where 


Di  is  the  thermal  diffusivity  of  the  gas  defined  as  Di  =  Ki/pjC  ,  Ki  is 

l 

its  specific  heat  at  constant  pressure.  The  polytropic  exponent  k  is 


given  by 


These  terms  are  also  contained  in  the  computer  program  used  to  evaluate 
the  theory. 


RESULTS  AND  DISCUSSION 


In  this  section,  the  major  aspect  of  this  work  is  examined.  The 
experimental  measurements  of  the  bubble's  position  in  the  sound  field, 
its  "rise-time",  and  the  voltage  measured  by  the  external  pill  trans¬ 
ducer  are  used  to  obtain  values  of  the  levitation  number  as  a  function 
of  a  radius  parameter  for  various  acoustic  pressure  amplitudes.  These 
values,  which  were  calculated  using  the  computer  program  in  Appendix  A, 
show  a  variation  of  the  levitation  number  over  a  range  of  normalized 
bubble  radii  for  both  linear  and  nonlinear  oscillations.  The  normal¬ 
ized  radius  is  the  value  of  the  bubble's  radius  divided  by  the  resonance 
radius  for  the  given  conditions. 

Many  data  sets  were  taken;  shown  here  are  some  representative 
examples  at  a  frequency  near  22  kHz.  All  data  presented  here  were  taken 
with  constant  values  of  the  acoustic  pressure  amplitude  and  the  tempera¬ 
ture  (for  a  given  data  set) ,  and  these  both  were  then  varied  over  a  wide 
range.  Some  representative  data  are  shown  in  Fig.  7  in  which  the 
temperature  was  35°  C  and  the  acoustic  pressure  amplitude  was  0.16  bars. 
The  symbols  are  experimental  points;  the  dashed  and  solid  lines  represent 
the  prediction  of  the  linear  and  nonlinear  theories  respectively. 

For  bubbles  less  than  40%  of  the  resonance  radius,  the  levitation 
number  is  nearly  a  constant.  These  points  agree  very  well  with  theory. 
Since  agreement  of  small  bubble  radii  with  theory  was  consistently  good, 
these  small  bubbles  were  used  in  each  data  set  to  calibrate  the  sound 
field  within  the  levitation  cell  in  terms  of  the  voltages  measured  by 
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the  pill  transducer.  This  section  of  the  graph  agrees  well  with  theory 
because  the  bubbles  are  just  oscillating  linearly  and  the  linear  theory 
coincides  with  the  nonlinear  theory  in  this  region.  In  the  area  from  40% 
to  55%  of  resonance  radius,  the  bubbles  clearly  go  into  nonlinear 
oscillations  as  indicated  by  the  rapid  increase,  then  decrease  of  the 
levitation  number.  The  acoustic  pressure  amplitude  for  this  data  is 
approximately  0.16  bars  and  the  levitation  number  varies  from  about  0.08 
for  small  bubbles  to  about  0.10  for  bubbles  near  harmonic  resonance. 
Bubbles  with  a  normalized  radius  greater  than  0.55  again  yield  levitation 
numbers  in  close  agreement  with  theory. 

Fig.  8  shows  data  using  the  glycerine/water  mixture  at  36°  C.  The 
acoustic  pressure  amplitude  was  0.19  bars  for  this  case.  This  higher 
pressure  indicates  the  relationship  between  the  levitation  number  and 
acoustic  pressure  amplitude  when  compared  with  Fig.  7.  As  the  acoustic 
pressure  amplitude  increased,  so  does  the  levitation  number.  The  match 
with  theory  for  small  and  large  bubble  radii  is  still  apparent,  as  well 
as  the  rapid  change  in  levitation  number  (up  to  about  0.11)  near  harmonic 
resonance. 

Fig.  9  shows  data  taken  at  39°  C  for  glycerine/water.  The  acoustic 
pressure  amplitude  was  0.17  bars  and  the  levitation  number  for  small 
bubbles  in  linear  oscillations  is  almost  constant  at  about  0.08.  This 
result  is  consistent  with  the  two  previous  figures.  The  large  bubbles 
with  a  normalized  radius  greater  than  0.55  also  lie  along  the  solid 
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Normalized  Bubble  Radius 

Variation  of  levitation  number  with  normalized  bubble  r 


Data  was  taken  for  the  glycerine/water  mixture  for  a  variety  oC 
temperatures,  some  as  high  as  62°  C.  Because  of  the  large  difference 
in  temperature  between  the  ambient  air  and  the  liquid  in  the  cell,  data 
was  very  difficult  to  obtain  at  the  higher  temperatures.  Convection 
currents  and  variable  gas  content  levels,  in  addition  to  evaporation, 
caused  many  delays  and  complications.  Shown  in  Figs.  10  and  11  are  data 
obtained  in  the  "high"  temperature  region. 

The  data  shown  in  Fig.  10  represent  a  bubble  in  glycerine/water  at 
49.5°  C.  The  acoustic  pressure  amplitude  was  0.16  bars.  The  levitation 
number  in  the  small  bubble  region  is  almost  constant  at  0.08,  indicating  as 
in  the  previous  data,  a  linear  dependence  on  acoustic  pressure  which  agrees 
with  the  theory.  The  levitation  number  in  the  range  of  harmonic  resonance 
changes  from  a  maximum  of  about  0.09  to  a  minimum  of  about  0.08,  which, 
when  compared  to  Fig.  7  with  the  same  acoustic  pressure  amplitude  of  0.16 
bars,  shows  a  considerable  reduction  in  range.  In  Fig.  7,  with  the  same 
pressure  amplitude  but  with  a  temperature  of  35°  C,  the  range  of  the 
levitation  number  over  the  harmonic  resonance  was  from  a  maximum  of  about 
0.10  to  a  minimum  of  about  0.08.  This  change  is  related  to  the  higher 
temperature,  which  among  other  things,  should  affect  the  damping.  Pre¬ 
liminary  calculations  indicate  that  the  effect  of  the  additional  vapor 
in  the  bubble  is  to  increase  the  damping.  For  bubbles  of  large  size, 
the  agreement  with  theory  is  not  as  good  as  in  the  data  taken  at  lower 
temperatures.  This  effect  may  be  due  to  the  difficulty  of  measuring 
the  correct  radius  of  the  larger  bubbles.  Large  bubbles  in  the  liquid 
at  relatively  high  temperatures  tended  to  grow  very  rapidly  to  resonance  size 
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and  so  data  acquisition  was  more  difficult  than  at  lower  temperatures. 

In  Fig.  11,  data  is  shown  for  a  temperature  of  50°  C.  In  this 
case  bubbles  of  large  radius  were  not  obtained;  however,  the  range  of 
the  levitation  number  in  the  region  of  0.4  to  0.5  normalized  radius  was 
from  a  maximum  of  about  0.14  to  a  minimum  of  about  0.10.  This  wide 
range  is  due  to  the  relative  large  acoustic  pressure  amplitude  of  0.22 
bars.  The  smaller  bubbles  in  this  particular  figure  show  a  wide 
distribution  of  levitation  numbers.  This  effect  was  probably  a  result 
of  the  currents  in  the  liquid  at  such  high  temperatures.  It  is  apparent 
from  the  data  presented  in  Figs.  7  through  11  that  the  levitation 
number  is  affected  considerably  by  the  temperature  of  the  host  liquid. 

In  addition  to  the  glycerine/water  mixture,  isopropyl  alcohol  was 
also  used  as  a  host  liquid.  The  frequency  used  was  near  20.5  kHz. 
Temperature  was  varied  from  9°  C  to  30°  C,  while  acoustic  pressure 
amplitudes  ranged  from  0.075  to  0.23  bars.  Fig.  12  shows  data  taken  a 
9°  C  with  an  acoustic  pressure  amplitude  of  0.08  bars.  As  with  the 
glycerine/water,  small  bubbles  in  the  range  of  0.3  to  0.4  for  normalized 
radius  were  used  to  calibrate  the  cell  for  each  data  set.  These  small 
bubbles  have  a  nearly  constant  levitation  number  of  0.04  which  agrees 
well  with  the  theory,  indicated  by  the  solid  curve.  In  the  region  of 
normalized  radius  from  0.4  to  0.55,  the  levitation  number  changes  from 
a  maximum  of  about  0.046  to  a  minimum  of  about  0..044.  These  data  in 
isopropyl  alcohol  are  similar  to  those  in  glycerine/water  in  that  again 
the  theory  predicts  a  much  smaller  variation  in  the  levitation  number  than 
is  observed,  even  for  this  case  in  which  the  amplitude  of  the  variation 
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evitation  number  with  normalized  bubble  radius  for  glycerine/water 
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is  small.  If  the  data  are  examined  away  from  the  harmonic  resonance 
region,  the  match  with  theory  is  quite  good.  At  this  low  temperature, 
it  was  possible  to  obtain  a  very  wide  range  of  bubble  radii.  The  two 
different  symbols  on  this  figure  refer  to  two  different  drag  laws 
that  were  used.  The  solid  circles  refer  to  a  solid  sphere  drag  law 
whereas  the  squares  refer  to  a  fluid  sphere  law. 

Fig.  13  shows  data  taken  for  a  temperature  of  15°  C  with  an  acoustic 
pressure  amplitude  of  0.12  bars.  The  linear  portion  of  the  graph  for 
small  bubble  radii  shows  a  near  constant  levitation  number  of  about 
0.056  which,  when  compared  to  the  other  figure  shown,  shows  that 
increasing  the  acoustic  pressure  amplitude  increases  the  levitation 
number.  The  experimental  data  again  agree  with  theory  for  small  bubbles 
and  large  ones,  but  is  not  in  agreement  in  the  region  of  the  harmonic 
resonance. 

In  Fig.  14,  the  data  was  taken  at  20°  C.  The  acoustic  pressure 
amplitude  was  0.12  bars  and  the  levitation  number  for  the  small  bubbles 
with  normalized  radius  less  than  0.4  was  a  nearly  constant  0.056  which 
agrees  with  the  levitation  number  for  the  same  region  in  Fig.  13  (same 
acoustic  pressure  amplitude) .  As  the  bubble  grew  to  the  region  of  the  harmonic 
resonance  the  levitation  number  again  changed  much  more  than  the  theory 
predicted.  After  the  bubble  reached  55%  of  resonance  radius,  there  was 
again  agreement  with  theory.  The  fact  that  the  magnitude  of  the  oscilla¬ 
tion  in  the  levitation  number  near  the  harmonic  resonance  is  considerably 
larger  at  20°  C  than  at  15°  C,  while  the  acoustic  pressure  amplitudes 
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remained  the  same,  indicates  that  the  presence  of  small  amounts  of 
additional  vapor  inside  the  bubble  significantly  affects  its  nonlinear 
response. 

Fig.  15  shows  data  taken  at  26°  C  with  an  acoustic  pressure 
amplitude  of  0.105  bars.  The  levitation  number  for  small  bubbles  in  the 
linear  region  is  near  0.05.  This  value  is  less  than  those  in  Figs.  13 
and  14  which  is  expected  due  to  the  reduced  acoustic  pressure  amplitude. 

It  is  difficult  to  compare  the  glycerine  and  water  data 

with  the  isopropyl  alcohol  data  since  the  acoustic  pressure  amplitudes 
were  not  the  same.  The  acoustic  pressure  amplitude  was  determined  in 
part  by  the  liquid  itself,  and  it  was  not  possible  to  select  a  value 
before  a  run  was  made.  Even  though  several  data  runs  were  made,  a  set 
of  data  investigating  the  variation  in  levitation  number  with  temperature 
for  fixed  acoustic  pressure  amplitude  was  never  obtained.  Furthermore, 
because  the  liquid's  physical  properties  have  an  effect  on  the  magnitude 
of  the  levitation  number,  comparisons  of  data  sets  for  the  same  values 
of  the  temperature  and  acoustic  pressure  cannot  realistically  be  made. 

All  figures  show  strong  agreement  between  experimental  data  and 
theory  for  the  regions  above  55%  and  below  40%  of  resonance  radius,  that 
is, in  the  linear  regions.  The  principal  area  of  concern  is  the  region 
near  the  harmonic  resonance.  In  Fig.  16,  a  typical  set  of  experimental 
results,  using  the  glycerine/water  mixture,  is  shown.  The  linear  theory 
prediction  is  shown  by  the  dashed  line  and  the  nonlinear  theory  showing 
a  combination  of  contributions  from  the  individual  resonance  regions. 
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is  indicated  by  the  solid  line  marked  with  Damping  x  1.0.  The  other 
solid  line  marked  Damping  x  0.25  shows  how  the  theoretical  values 
changed  when  the  damping  was  reduced  by  a  factor  of  4.  This  results 
in  a  very  close  fit  with  experimental  data.  The  areas  which  agreed  before 
continue  to  do  so  and  the  region  of  principal  concern  now  matches  nicely 
with  theory.  The  fact  that  this  simple  adjustment  in  the  damping  creates 
such  good  agreement  indicates  a  flaw  in  the  theory.  Although  it  would 
be  easy  to  attribute  this  discrepancy  to  an  incorrect  evaluation  of  the 
damping,  such  an  assumption  is  probably  unwarranted.  Independent 
measurements  of  bubble  damping  show  no  discrepancy  with  theory,  and  to 
expect  one  here  is  not  reasonable.  What  is  probably  occurring  is  that 
the  polytropic  exponent,  which  is  treated  with  a  linear  theory,  must 
also  be  examined  nonlinearly,  in  order  to  assign  it  an  appropriate  value 
in  this  highly  nonlinear  region.  Also  shown  in  the  figure  are  the 
variations  of  the  phase  angle  ^  as  a  function  of  radius .  The  rapid 
change  in  the  phase  angle  in  the  region  of  harmonic  resonance  is  the 
main  cause  for  the  rapid  change  of ^evitation  number  near  the  normalized 
radius  value  of  0,5. 

It  is  clear  that  work  needs  to  be  done  in  the  theoretical  area. 

The  experimental  measurements  required  to  evaluate  the  levitation 
number  are  quite  straightforward,  and  thus  there  is  little  cause  to 
expect  the  values  of  the  levitation  number  to  be  in  error.  An  examina¬ 
tion  of  the  theory  in  a  completely  nonlinear  way  is,  however,  a  very 
difficult  task  and  cannot  be  done  quickly.  Yet,  these  data  show  that 
some  interesting  and  important  effects  are  being  observed  that  offer 
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Errors  in  this  experiment  could  rise  from  several  areas.  In  the 
direct  measurements  involved  in  each  data  set  —  the  temperature, 
frequency,  bubble  position,  rise  velocity  and  acoustic  pressure  —  are 
of  major  concern.  In  addition  to  these  values,  physical  constants 
associated  with  the  liquids  and  gases  were  utilized.  In  some  cases, 
these  values  were  obtained  from  references  and  in  other  instances,  the 
values  were  measured  experimentally.  Examinations  of  each  area  for 
possible  error  follow. 

The  values  used  as  constants  were,  in  general,  very  accurate  for 
the  temperature  used.  The  difficulty  arose  in  the  fact  that  temperature 
was  measured  to  an  accuracy  of  0.2°  C,  giving  an  error  of  0.7%  over  the 
range  of  temperatures  used.  If  the  viscosity  for  isopropyl  alcohol  at 
30°  C  was  obtained  from  a  table,  this  value  was  valid  for  that  temperature 
but  the  temperature  of  the  liquid  could  have  been  different  from  measured 
values  or  could  have  changed  over  the  period  of  data  acquisition.  To 
avoid  error  related  to  this  circumstance,  temperature  measurements  were 
made  before  and  after  data  collection  and  data  was  taken  over  as  short 
a  time  as  possible,  on  the  average  of  30  minutes.  In  order  to  deter¬ 
mine  the  el  feet  on  the  levitation  number  or  the  normalized  radius  of 
an  inprecision  in  one  of  the  many  constants  and  variables  used  in  their 
calculation,  a  sensitivity  analysis  was  performed.  That  is,  the  parameter 
under  consideration  was  varied  by  an  amount  equal  to  its  inprecision  and 
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then  the  applicable  computer  program  was  utilized  to  test  the  effect 
of  this  variation  on  the  levitation  number  or  normalized  radius.  For 
example,  viscosity  was  measured  to  an  accuracy  of  ±0.1%.  Substituting 
the  two  values  of  viscosity  into  the  computer  program  resulted  in  a 
change  in  the  normalized  radius  of  0.04%  and  had  very  little  effect  on 
the  levitation  number.  Thus,  the  levitation  number  and  the  normalized 
radius  were  rather  insensitive  to  errors  in  the  measurement  of  the 
viscosity. 

The  values  relating  to  the  properties  of  the  acoustic  wave  in  the 
cell  were  determined  by  calibration  of  the  cell  with  a  probe  hydrophone. 
The  errors  involved  in  these  values  were  a  result  of  various  factors. 

The  probe  itself  tended  to  change  the  sound  field  and  the  level  of  the 
liquid.  Other  considerations  related  to  the  voltage  and  position 
measurement  will  be  mentioned  later.  These  factors  were  counteracted, 
to  a  degree,  by  the  use  of  a  calibration  routine  using  small  bubbles 
in  linear  oscillations.  This  was  done  for  each  data  set.  The  measured 
wavelength  of  the  acoustic  signal  was  a  critical  factor  and  was  accurate 
to  only  1.1%.  This  results  in  a  variation  of  the  levitation  number  by 
a  value  of  1%. 

Another  area  of  concern  relating  to  the  error  in  values  used  for 
constants  is  that  when  the  mixture  of  glycerine  and  water  was  used, 
values  of  the  viscosity  were  obtained  from  measurements  already  mentioned 
These  measurements,  even  though  they  fit  well  in  between , viscosity  values 
of  known  mixtures  were  probably  less  accurate  than  textbook  values. 
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In  the  area  of  measured  values,  the  frequency,  which  was  continually 
tuned  to  the  correct  resonance  value  ,  was  accurate  to  within  0.05%. 

This  error  resulted  in  no  measurable  effect  on  the  calculated  radius  or 
the  levitation  number.  The  values  of  rise  time  were  accurate  to  0.01 
seconds  as  a  result  of  the  methods  used  to  obtained  measurement.  The 
bubble  was  positioned  in  the  field  of  view  and  its  position  was  recorded. 
Then  the  sound  was  discontinued  and  the  bubble  was  allowed  to  rise  a 
certain  distance,  as  indicated  by  the  fiducial  lines  in  the  microscope. 

It  could  be  determined  that  the  bubble  was  not  changing  size  by  the  fact 
that  rise-times  could  be  repeated  within  0.02  seconds  several  times. 

This  inaccuracy  of  the  time  resulted  in  a  maximum  0.3%  error  in  the 
calculated  normalized  radius.  An  effort  was  made  to  keep  rise  times 
around  2  seconds  since  longer  times  could  result  in  a  change  in  bubble 
radius  during  its  rise  while  shorter  times  increased  the  relative  error. 

The  bubble  position  was  determined  by  use  of  the  cathetometer 
(accurate  to  0.02%).  This  reading,  when  related  to  absolute  displacements 
within  the  levitation  cell,  would  give  information  relating  the  position 
of  the  bubble  in  the  sound  field  to  the  distance  the  bubble  traveled 
during  the  given  rise  time.  The  location  of  the  bubble  in  the  sound 
field  was  used  to  determine  the  bubble  radius  to  an  accuracy  of  0.03% 
and  this  information  was  used  to  give  a  levitation  number,  accurate  to 
0.07%. 

Bubble  displacement  as  measured  through  the  use  of  fiducial  lines 
in  the  microscope  was  accurate  to  0.1  line  which  gives  an  accuracy  of 
1.2%  since  the  bubble  was  allowed  to  rise  an  average  of  2  seconds.  This 
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translates  to  bubble  movement  of  between  1  mm  and  16  mm.  The  accuracy 
results  in  a  radius  variation  of  1.4%. 

Voltage,  which  was  a  measure  of  the  acoustic  pressure  amplitude  in 
the  cell,  was  held  constant  during  data  collection  but  was  periodically 
checked.  The  voltage  would  only  change  when  the  bubble  grew  to  resonance. 
An  accuracy  of  0.001  volts  results  in  0.2%  error.  The  fluctuation  of 
voltage  was  an  indication  of  the  presence  of  additional  bubbles.  These 
bubbles  usually  were  levitated  near  another  antinode  or  on  the  surface 
of  the  cell.  If  these  bubbles  were  present,  they  would  tend  to  disrupt 
the  sound  field  if  they  grew  to  resonance  size.  Data  was  not  taken  when 
the  additional  bubbles  were  present. 

The -calculations  of  bubble  radii-and- levitation  number  are  very 
complex  as  shown  by  the  computer  program  1  in  Appendix  A.  The  resulting 
errors  in  experimental  data  could  be  as  high  as  20%  in  some  instances, 
but  the  values  shown  are  in  such  strong  agreement  with  theory,  in  the 
applicable  regions,  that  it  is  felt  that  a  5%  error  would  be  more 
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CONCLUSIONS 


A  unique  and  sensitive  technique  has  been  developed  for  the  study 
of  nonlinear  oscillations  of  an  individual  bubble  in  liquid.  Measure¬ 
ment  of  various  aspects  of  the  individual  bubble's  behavior  is  possible 
with  this  levitation  technique.  The  technique  has  been  very  successful 
in  the  limited  study  explained  in  this  paper  and  much  promise  is  shown 
for  future  study.  These  measurements  may  be  made  with  little  difficulty 
after  some  experience ' is  gained  using  this  technique. 

The  fit  of  experimental  data  with  theory  is  consistent  throughout 
the  various  data  samples  shown  for  the  regions  of  normalized  radius 
less  than  0.4  and  more  than  0.55.  Theory  and  experiment  do  not  agree 
in  the  region  of  the  n=2  harmonic  unless  the  damping  is  reduced  by 
a  factor  of  4.  The  reasons  for  this  discrepancy  are  unclear  and  more 
study  is  needed  in  this  area. 
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APPENDIX  A 


COMPUTER  PROGRAM  1 
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INTEGER  N>  CONST  .  J5>  05 

REAL  A8.B8.C8.E8.F8.G8.R9.  B9.  C9.  D9.  ETREX.  SLO 

RERL  Q6 .  W4 > M5  >  W6  .  V6 . td?  »  V7.  W8 .  V8 » R9 

REAL  EE » EP  »  YE .  YP » PE  C2  0 0>  .  PR  <200>  .  WE  <20G>  >  UP  <200> 

RERL  PE3  C200)  .  PP3  <200)  »  C6  <200>  >  CON.  PRE  C200>  >  RP  <200> 

PERL  B1 . B2 >  B3.  S5  <200>  .E5<200>  .  P5  <200>  »  R7  C200)  , R2  <200) 
RERL  A5»  F5  <.'200>  .H5<200>  »G5<200>  >  K. R6. R8  <200)  » WO*  RRT »  ETR1 
RERL  F7< 20 0> 

RERL  SUR. DEL. PZ. ET.  ER»  DIF?  ETR2 
RERL  C0«  Ml  .  M2.L5. SI  .  T1  *  Q1 ,  T5  <200>  .  C5  <200>  .  V  <200>  .  D5  (200> 
PEAL*8  M»  D1 .  T.  F«  G.  R1 .  W>  G1  •  X 

REALMS  0.  P0» Cl .  DV.  DR.  DT -XI .  DTI  *  DT2.  El 
REALMS  BT3.  E2»  BE.  DP?  ETA*  C2 
REALMS  G2.  A«  B»C«  D.  E»  H.  I  .  J.  L.  VI *  Y2>  Y 
REALMS  Z1  .  Z.  Z2.  Q.  fc  AP.  DPV.  DPR.  DPT 
FNC  OO  =  CDEXP  CO  +  BEXP  <-X>  ->  /2 
FNS  <X >  =  <  DEXP  <X>  -DEXP <-X> > /2 
PERD  <9.  1  0>  CON.  Ml  >  M2. L5 
RERD  C9. 13>  Sl.Tl.Q5 
PERD  <9 .  I  0>  B 1 .  B3 »  B5 «  U 
PERD  <9 .  1  0-'  G .  F .  C 1  .  C2 
PERD  <  9 *  1 0>  M . D 1 . R9 . SUR 
250  FORMAT  <  •'  123.  LAM.  «.  CO  AND  SLO'> 

220  FORMAT <I> 

390  FORMAT C2F> 

400  FORMAT </  L5  AND  CO  =  > 

1 0  FORMAT  <4F> 

13  FORMAT < £F. 1I> 

17  FORMAT <3F> 

18  FORMAT  '2F> 

TYPE  250 

ACCEPT  220. CONST 

ACCEPT  300. L 5 

ACCEPT  220. Q5 

ACCEPT  300.  CO 

ACCEPT  300. SLO 

IF  'CONST. LE. 2>  GO  TO  33 

TYPE  312 

312  FORMAT ('  N  RADIUS  ETREX 

1  ETR-ITE  KAPPA' > 
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83  CONTINUE 

IF  CCONST.  GE.  3>  GO  TO  22 
TYPE  310 

310  FORMAT  (.'  N  NOR-RAD  P-ELLER 

IALPHA-PRI'> 

22  READ  <9*  120)  <T5  <N>  ?  C6  <‘N>  »  V  <N>  t  D5  <N>  » N=1 » Q5> 

120  FORMAT <4F> 

300  FORMAT  <1F> 

P 0=1. 0E6 
T=2?3+T1 
Rl=8. 3D7 

190  K=6.23 32/L5 

U=2. *3. 14159*F 
G1  =M^D1  ♦W-"'  <G4G^R1*T> 

DO  130  N=1*Q5 

C5  <:n>  *  <C6  cn>  -c  o>  ♦m  1 

S3  <:n>  =D5  <:n> 

E5<N>=S5<N>/T5CN> 

R5<N>=1. OE-3 

36  R?  <N> =2. »R5  CN> »E5  <N> ♦CE^'U 

R2  >'N>  =1 .  +0.  19?*R7  <N>  ♦♦O.  63+ <2.  6£-4>  ♦R?  <N>  *^1 . 38 
R5  <N> =SGRT  <4. 5*U*E5  CN> +R2 <N> ✓ < 930 ♦C£> > 

A5=SQRT  <  <R5 <N> -R6> ##2> 

IF  CA5.LE. 1. 0E-5>  GO  TO  40 
R6=P5  <N> 

GO  TO  36 

40  CONTINUE 

R3<N>«1. 000*R3(N> 

F3  <N>  =B1  ♦SIN  CIOC5  <N>  >  +B3»SIN  <3  ♦IOC  5  <  N>  > 

1  +B5^SIN<5.  ♦K>C5CN>> 

H5  01>  =3^B5^C0S  <S^k>C5  <N>  > 

G5  <N>  =K>  <B1 ♦COS  <K>C5  <N>  >  ♦3^B3^C0S  <3  ♦IOCS  <N>  >  +H5  <N>  > 
G2*W*R3 <N> ♦RS •N> sDtsG 
A=  * <G-1  >  «'G>  ♦  <Gl  'y2.  > 

B=DSQRT  <GOG2> 

C=DSQRT  <  G*G2-'  2.  > 

D=2. ♦€♦ <1-A> 

E=2. ♦€♦ <1+R> 

H=FNC  <  D> -DCOS  <  E> 

I =FNC  <2. ♦A^B> -DCOS  <2. ♦B) 
j«B#<A#G1/G-1> 

L=B»<A+G1'G> 

Yl=2. ♦A^I  ♦ <D*FNS <  D>  +E»DSIN<E) -2. »H> 

Y2=2. ♦fOCJ^FNS  <2. ♦A^B> -MODS  IN  <2. +B> -l+Gl'G> 

Y=Y1-Y 2 

Zl=2. ♦A^O<E»FNS  CD)-D*DSIN<E)  > 


Z2=2. +H+ CL*FNS C2.  -J»DSINC2. *B> -I> 

Z=21-Z2 

Q1 =2. ♦ <G1 *Ys 6-2. ♦R^Y +Z> 

Q=Q1>"CY*Y+Z^Z> 

KRP=G^G1 ♦G2^Q^3. EO 

X1=R5  CN)  ♦DSQRT  <2.  ♦W^DO 

DT1*X1*<FNS  CXI)  +DSIMCXD)  -2.  ♦  CFNC  CXI) -BCDS  CXI)  ) 
DT3=FNS CXI) -DSIN  CXI) 

DT2=X1  #X  1  ♦  CFNC  CXI  >  -DCOS  <X  1  )  >  +3.  ♦  CG- 1 .  >  ♦X 1  *013 
DT=3. ♦  CG-l . > ♦DT1/DT2 
E2=  CFNC  CX 1 > -DCOS CX 1 >  >  #X 1 

El =3. ♦ CG-l) ♦ CFNS  CXI) -DSIN CXI) > /E2 
ETR=G'  <  <1+DT+DT>  ♦  Cl.+El)) 

BR=C£»R5CN>  ♦♦S.  ♦U*^3.s<3.  ♦ETR^P0*C1> 

DV=4.  ♦W^U''  C3. ♦ETfl*PO> 

BE=DV+DT+BR 

DP=  CY-Z^Gl'G-2.  ♦R^Z)  ✓  (Gl+Y/6-2.  ♦H*Y+Z> 

DPV=DV 

DPR*DR^  < 1 . ♦ <W*R3 CM) ✓ C 1 >  ^*2 . > 

DPT=DPV+DPR+BP 
EE=SGRT CC2/3. /ETfl/PO> *R5  <N> 

EPsSQRT  CC2>'3.  ✓  KRP'PO)  ♦RS  CN)  ♦« 

YE*  C  <1  -EEp*2. >  ♦♦2.  +DE^2>  ✓  Cl  -EE*^2.  > 

YP*  C  Cl -EP»#2.  >  ♦♦£.  +BPT^2.>/C1-EPpp2.> 

PE  cn>  =SQRT  C2^ETfl^C2*98i>P0.'F5  CN)  sG5  CM)  ♦YE> 

PP CN)  *SQRT  C2*KRP^C2*980PP0^F5  CN)  '65  <N>  ♦YP) 

PE  CN)  *1 . 0E-6*PE CN) 

PPCN)*1. 0E-6*PPCN> 

WE  CN)  =V  CM)  .'PE  CM)  /SLO 
MP  CM)  =V  CM)  sPP  'N)  ^'SLD 
PRE  CM)  *  <  v  CN)  >'SLD)  »P0 
fiP  CM) *£♦€£♦980''  <3*PRE CN) *G5  CM)  > 

PE3cN)  *PE  *N)  ♦FS  <  N) 

PP3 CM) =PP  CN) »F5  CM) 

F? CN) *RP  CM) '  PE3  CM) 

R8*C2*M*W*R5 CM) *R5  CN) / C3*P0) 

B8»BT 

C8*C£p CR5  CN) ♦♦S) ♦<W*^3) ✓ C3*P0»C1) 

Bo=4*M»L)/  C3»P0> 

ESs2+C2+980.  *POs  CF5  CN)  *G5  CN)  ) 

G8=C8+D8 

P8*CVCN)xSL0>  ♦PO 

F8* CP8*P8) ✓ E3 

R9«1.+B8*B8 

B9*£*68*B8-£#R8-F8 

C9=G8^G8+R8»R8+F8*R8 

D9»SQRT  CB9*B9'  <4*R9*R9> -C9^R9) 
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R8  (N>  =R 5  (N> ✓  R9 
PZ=P0+2+SUR' R5 (N> 

ET=2+SUR' (3+PZ*R5  (N>  > 

ETR2=KRP 

611  WO=SQRT  < (3+ETR2+PZ-2+SUR'R5  (N>  > ✓ (C2+R5  (N> ♦♦a)  > 
BEL=RTRN(W»WO+BPT/'  (WO+WO-W+t4>  > 

rrt=  a .  -«♦«✓  cwo+wo>  >  ♦♦2+  (w+bpt^uo>  ♦♦a 

ER=2+C2+980+PZ«>SQRT  (RRT) 

ETfll=P8*+2*F5(N>  *65  (N>  +COS  CBEL>  ' Efl+ET 
BIF=SGRT((ETR2-ETRl>++2> 

IF  (DIF. LE. 1 . 0£-3>  60  TO  612 
ETR2=ETfll 

612  ETREX=D9-B9^ (£^R9> 

IF  <C0NST-2>  230*  ISO*  290 
290  TYPE  130*N*R8CN> * ETREX* ETR1 ,  KRP 
IF  (CONST. EQ. 3>  GO  TO  180 
280  TYPE  1 30* N* R8  CN> *  PE 3 CM>  >  RP CN> *  V  (N> *  ME CN3 
130  FQPMRT (1 1  *  5F1  0. 6) 

180  Q6=Q5 

W4=  0 
W3*0 
V4=0 
V5=0 
V7=0 
M  7=0 

IF  ( CONST. EQ.3>  GO  TO  350 
BO  140  N=1*G6 
U4=W4+ME (N> 

140  V4*V4+WP (N> 

1,16=1.14^136 
V6*V4^G6 
BO  150  N=  l*i36 
1.17=1.17+  (ME  (N>  -W6-)  ♦♦£ 

150  V7=V7+  (WP  < N>  -V6>  ♦♦£ 

WS-SQRT  <UI7/’  (  Q6- 1  >  > 

V8=SGRT  ('■/?/  (Q6~l>  > 

TYPE  1 70*  CO*  W6*  W8 
IF  (CONST. EQ.2>  GO  TO  420 
TYPE  400 
TYPE  390*L5*C0 
170  FORNRT  <3F1 4.  6> 

200  FORMRT  <1 1 » 1F12.6) 

420  C0=C0-1 

IF  (CONST. £0. 1 >  GO  TO  350 
IF  (L5. LE. 29. >  GO  TO  130 
IF  (CO.  LE.  0.  >  GO  TO  190 
350  CONTINUE 
STOP 
END 
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THIS  PROGRAM  COMPARES  OUR  EXPERIMENTAL  DATA  WITH  THE  THEO 
RY 

OF  ANDREA  PROSPERETTI.  IT  CALCULATES  A  BUBBLE'S  MAXIMUM 
CHANGE  IN  RADIUS  AS  A  FUNCTION  OF  FREQUENCY  FOR  RESONANCE 
AND  THE  FIRST  AND  SECOND  HARMONICS. 


COMMON  AA2 

COMMON/ WR/AMPL <3>  » COSPH 
REAL  XMAX1 <3>  » XMAX2  <3> 

DATA  W 0 « W 1 >  W2 • XM AX 0 *  XMAX 1  *  XMAX2/1 0^0. / 
RRES= 0.0174 
P=0. £0* RRES 

♦♦♦♦  INITIALIZING  CONSTANTS 
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GAMMAP=1 .40 
BIGD1=. £05 
SIGMA=2£. 0 
PI NF=1 000000. 0 
RQE=. 7800 
DELTR0=. 01 
DELTR1=. 005 
TYPE  137 

FORMAT <'  INPUT  PA') 

ACCEPT  135?  PA 
FORMAT  (1F> 

TYPE  136 
FORMAT <.'  R/RO 
DFUDGE=0.  000 

BIGO=20£00^2»3. 141596 
SPEED* 161 000. 0 
AMU=  0. 022 
DO  100  1=1 .300 
PROSP=R/PRES 
PRQSP=R/RRES 
IF  <PROSP. GT. 0. 76)  GOTO  50 
ECKS=R^SQPT  <2. ♦BIG0/BIGD1> 

QQZ1  =  <SINH  CECKS) -SIN  <ECKS>  > / <COSH <ECKS) -COS  <ECKS>  > 


LEV  NUM 


ALPHA  PRIME') 


QQZ2=  <SINH <ECKS>  +SIN  (ECKS>  > ✓ <CDSH  <ECKS> -CDS  (ECKS)  > 
DTHERM=3*  (GAMMAP-l > ♦ <ECKS*QQ22-2>  ^ECKS 
BTHERM=DTH£RM/ (ECKS+3*  (GAMMAP- 1 > ♦QQZl > 

GAMMA=GAMMAP/  <  1  +DTHERM^2> 

GAMMA=GAMMA>'  <  1  +3*  <G AMM AP- 1  >  ♦QQZl /ECKS> 

6AMMA=1 .  O^GAMMA 
PO=PINF+2.  ♦SIGMfl/R 
BBLU=2. ♦SIGMA;' <R^PO> 

ALPHR1=4.  S^GAMMA^  (GAMMR+1 . > -2. ♦DELU 

ALPHA2=.5^GAMMA^<9.  ♦Gflf1f1fl^*2.  -*-18.  ♦GAMMA+l 1 .  >  -3.  ♦BBLU 
BIGDO=SQRT  <  <3 .  ♦GAMMA^P 0-2 .  ♦S I  GMR/R>  ✓  <ROE^R^R>  ) 
DMEGAO=R^BIGOO^SQRT<ROE/'PO> 

OMEGA=R^BIGO^SQRT  (ROE/PO> 

DVIS=4. ♦BIGD^RMU/ <3. ♦GAMMA^PINF) 

DRAD=R0E^ (R^B I GD) ♦♦ 3 ✓ <3 . ♦6RMMfl*P I NF^SFEED> 

DTDTAL=DTHERM+DVIS+BRAH+DFUDGE 

BCHECK=PTOTAL-DFUDGE 

B=. 5^R^BI6D^DTDTRL*SGRT  (ROE'  P0> 

B=B*BIGOO»BIGOO/ <BIGD^BIGQ> 

ETA=PA^PIMF 

XI  =  <1 .  -DBLU>  ♦ETf. 

RAZ=  <  (DMEGA^2. -DMEGA0^2. > ♦♦£.  +  C2. ♦B^DME6A> ♦♦2.  > ♦♦-.5 
WRITE  (23*20>  DTHERM..  DV1S,  DRRDt  DFUDGE*  BTDTAL*  DCHECK*  B 
20  FDRflflT  (.'  '  *  3X»  7  <E8. 3»  3X)  > 

♦♦♦♦  CALLING  SUBROUTINES  HRMNC1 *  HRMNC2»  RNB  RESNCE  ♦♦♦♦♦♦♦♦♦♦ 

WRITE  <22*  1 1 0>  _R»  PROSP.»  B>  DBLUji  ETR»  XI «  OMEGAO 
10  FDRMAT  <.'  '  *  E8.  3?  6  <2X»  E8.  3>  > 

IF  (PROSP.GE.  .20.  AND.  PROSP.  LE.  .320>  GD  TD  200 
60  IF  (PROSP.GE.  .  320.  RNB.  PROSP.  LE.  .600>  GD  TD  300 
70  IF  (PROSP. GE. . 600. RND. PROSP. LE. 1. 001)  GD  TD  400 
GOTO  140 

:00  CALL  HRMNC2  (ALPHA  1 >  ALPHA2.-  B*  OMEGAO..  XI ..  W8 •  XMRX2 •  IERR£»  ONE 
GA> 

DELTAR=DELTR0 
XMRX2 (3) =1 . 5^XMAX2  <3> 

AMAX2=XMAX2  <3> ✓  1 . 5 
TYPE  41 0» PROSP *XMRX2 (3>  » AM AX 2 
10  FDRMAT <'  ,»3X»E10.5»4X:»E10.5«4X-»3<E10.5» 3X> > 

GOTO  160 

00  CALL  HRMNC1  (PHASE.*  PHI » ALPHA1  >  ALPHA2»  B>  OMEGAO*  XI » W1  * 

/XMAX 1 » I ERR 1  *  DMEG A> 

DELTAR=BELTR1 

XMAX1 <3>  =1 . 5^XMAX1 <3> 

AMAX 1 SXMAX 1 <3>/1.5 

TYPE  41 0» PROSP > XMAX1 <3> » AMAX 1 
GOTO  170 


400  CALL  RESNCE  <AL  JHA1  *  ALPHAS »  B»  OMEGAO*  XI  .*  WO*  XMAX0»  OMEGA) 
DELTAR=DELTR0 
XMAX0=1 .S^XMAXO 
AMAX0=XMAX0/1 .5 
TYPE  41 0*  PRDSP*  XMAXO*  AMAXO 
410  FORMAT ('  ' » IX* 5  (El 1 . 5* 3X) ) 


IF  (IERR2.  NE.  0)  WRITE  <21 » 120)  IERR2 
IF  (1ERR1.NE.  0>  WRITE  <21*1 30)  IERR1 
120  FORMAT  <'  ♦♦ERROR  NOTE  IN  HRMNC2^*  '*I3> 

130  FORMAT  ('  ♦♦ERROR  NOTE  IN  HRMNC1^+  ',I3> 

40  R-DELTRR*RRES*R 
100  CONTINUE 
50  STOP 
END 

•♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦ 

THIS  IS  A  SUBROUTINE  TO  BE  USED  WITH  THE  PROGRAM  COMPARING 
OUR  RESULTS  WITH  THOSE  OF  PROSPERETTI .  IT  CALCULATES  THE 
RESONANCE  CURVE  FOR  THE  SECOND  HARMONIC  REGION. 

♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦ 

SUBROUTINE  HRMNC2  <ALPHA1 *  ALPHA2*  B*  QMEGAC*XI  *  Y»  XMAX*  I  ERROR* 
OMEGA) 

REAL  XMAX <3>  » XOFTAU  <2400) *  COEFFC  <4>  *  ROOT SR  < 3> * 

✓  ROOTS I (3) 

COMPLEX  ROOTS <3> 

COMMON  AA2 
NDEG=3 
I ERROR* 0 


Y=OMEGA^OMEGAO 

D=1 .  (OMEGA  0^2 .  -0MEGA^2.  > 

C0*(D^ (ALPHA 1-4. 5^0MEGA^2.  >  5)  ✓  (OMEGA  0^2 .  -4.  ♦0MEGA^2 

.  > 

Cl*  <.  5*D^ (D^  <ALPHA1  -1 . 5^0MEGA^2.  >  -1 .  >  >  /OMEGAO~2. 

C2*.  5^  <ALPHA1  +  1 3.  5^0MEGA^2.  >  ^  <DMEGA0^2.  -36.  ♦QMEGA^2. ) 
C3*<.  5»D^ <©♦  <ALPHA1  +1 . 5»0MEGA^2.  >  -1 . >  >  ✓  <DMEGA0^2.  - 
/  4. ♦0MEGA^^2. > 

C4*  <D^< ALPHA 1+4.  5+0MEGA++2.  >  -.  5>  ✓  (OMEGA 0++2.  -16.  ♦OMEGAS 


o 


l 
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C5=  (.  5^  (ALPHA1-13. 5^0MEGA^E.  >  >  /OMEGA  0^2. 

60=RLPHA1^(2.  ♦C5+C2)  -*-13.  5^0MEGA^2.  ♦  (.  25-2.  ^C2>  Z5*ALPH 

A2  ! 

Gl*-2. ♦ALPHA 1 ♦Cl +C4^ (. 5-D*  CALPHA1 -6. ♦QME6fl^^2. >  >  + 

/  1.5^(ALPHA2-.  5*aMEGR^*2.>  ^♦♦2.  -D 

G2=C0»  (D^ (ALPHA1 +3.  ♦0MEGR^^2. > 5> 

G5=C3^(D^(ALPHA1  +3.  ♦0MEGfl^^2.  >  -.  5>  +.  25^  (1 .  -D^  (ALPHA2+ 

/  1 . 5*DMEGA^2.  )  >  ♦D^^2. 

QUE3=9.  ♦OMEGA^E.  -QMEGA0^2.  -XI ♦♦E.  ♦  <G1 -G2) 

♦♦♦♦♦  INITIALIZING  THE  COEFFICIENTS  FOR  THE  CUBIC  "C^2. "  ♦♦♦♦♦ 

COEFFC  (1>=G0~2. 

COEFFC  (2) =2. ♦QUE3^G0 

COEFFC  (3)  =QUE3~2.  +  <6.  ♦OMEGR^B)  ♦♦E. 

COEFFC C4> =-l .  +G5++2. ♦XI ♦♦G. 

♦♦♦♦♦  FINDING  THE  ROOTS  OF  "C^»2. ”  ♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦ 

CALL  ZPOLR  (COEFFC*NDEG,ROOTS, IER) 

IF  (IER.NE.O)  IERROR=IER 

♦♦♦♦♦  LOOP  TO  FIND  THE  REAL  ROOTS  OF  C^2  AND  COMPUTE  THE  MAX. 
♦♦♦♦♦  VALUE  .OF  XOFTAU  FOR  EACH  OMEGA. 

DO  99  J=1 t 3 

ROOTS  I  (  J>  =  A I  MAG  (ROOTS  (  J>  > 

ROOTSR (  J)  =RERL (ROOTS ( J> > 

IF  (ROOTSI (J) . EQ. 0. >  GO  TO  98 
GO  TO  99 
98  KSTOP=J*80G 

KSTART=KSTOP 

X 1 =-6 . ♦QMEGA^B^SQRT (ROOTSR (J)  > ✓ <G5«XI«^3.  > 

X2=-SQRT (ROOTSR (J> > ♦ (OUES+GO^ROOTSR (J> > / (G5^XI~3.  > 

PHI =ATAN2 (X 1  *  X2> 

TAU=0. 

R=.  04.'  OMEGA 

DELTA=ATAN  (  (2.  ♦OMEGA^B)  ✓  C0MEGR~2.  -OMEGAO~2.  >  > 

♦♦♦♦♦  LOOP  THAT  CALCULATES  XOFTAU.  ♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦ 

DO  88  K=KSTART » KSTOP 
THETA=3. ♦OMEGAS TAU+PHI 

XOFTAU  (K) =SQRT (ROOTSR ( J> > ♦COS (THETA) +XI ♦ < (0MEGA^2. - 
OMEGA 0 ! 

/  ♦♦E. > ♦♦E. ♦ (2. ♦B^OMEGA) ♦♦E. > ♦♦-. 5^C0S (OMEGA^TAU+DELT 

A>  +  (C1  + 
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✓  C3^C0S  <2.  ♦OMEGA^TAU)  )  ♦XI^E.  +  <C5+C2^CDS  (2. ♦THETA)  )  ♦ 
ROOTSRj 

✓  <  J)  +  (C4^C0S  (OMEGA^TAU+THETA)  +CO*COS  (OMEGA^TAU-THETA 

>  >  ♦xi*J 

/  SORT  (ROOTSR  <J>  > 

TAU=TAU+R 
88  CONTINUE 

♦♦♦♦  LOOP  THAT  FINDS  THE  MAX.  VALUE  OF  XOFTAU.  ♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦ 

XMAX (J)  =XOFTAU  (KSTART) 

DO  77  L=KSTART,kSTDP 

IF  (XMAX(J)  .LT.  XOFTAU  <L))  XMAX (J) =XOFTAU (L) 

77  CONTINUE 

XMAX (J) =XMAX  (J) ♦COS  (DELTA) 

BJERK=AAZ^XI^COS (DELTA) 

XMAX  (  J>  =BJERK 
99  CONTINUE 


RETURN 

END 


THIS  IS  A  SUBROUTINE  TO  BE  USED  WITH  THE  PROGRAM  COMPARING 
OUR  DATA  WITH  THE  RESULTS  OF  PROSPERETTI .  IT  CALCULATES  T 


RESONANCE  CURVE  FOR  THE  FIRST  HARMONIC  REGION. 


SUBROUTINE  HRMNC1  (PHASE, PHI , ALPHA! , ALPHA2, B,  OMEGAO, XI , Y, X 
MAX, 

/I ERROR, OMEGA) 

REAL  XMAX  (3)  ,  XOFTAU  (2400)  ,  COEFFC  (4)  ,  RODTSR  <  3'.  , 

''  ROOTS  I  (3) 

COMPLEX  ROOTS (3) 

COMMON-'  WR/CC  (3)  ,  X2 
COMMON  AAZ 
NDEG=3 


I  ERROR= 0 


Y=OMEGA«'OMEGAO 


BIGD=1 . ✓' (OMEGA0^2. -OMEGA^2.) 

C0= (BI6D^(ALPHAl-3. ♦0MEGA#^2. ) -. 5) / CDMEGA0^2. -OMEGA^E. 

Cl* (. 5^BIGD^ (BIGD^ (ALPHA1 -1 . S^OMEGA^E. )-!.)) ^OMEGA0^^2. 
C2*  (.  $♦  (ALPHA1+6.  ♦OMEGA^E.  )  )  /•  COMEGAO^E.  -16.  ♦0MEGA^2.  ) 
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C4=  (BI GD+  CALPHA1  +3.  +0MEGA++2.  >  5)  ✓  (OME6AO++2.  -9.  +QME6A+ 

♦2.  > 

C5=  (. 5+<ALPHAl-6. +0MEGA++2.  ) > ^ OMEGA 0++2. 
C6=2.+0MEGA+BIGD++2. 

BETA2=.  5+BIGD+(BIGD+(ALPHAl+l .  5+  (0MEGA++2.  >  )  -1 .  > 
GZER0=ALPHA1 ♦  <2.  ♦C5+C2)  +6.  +0ME6A++2. +(.  25-2.  +C2>  75+ALP 
HR2 

G0NE=-2.  ♦PtLPHR  1 +C 1 +C4+  <.  5-BIGD+  (RLPHA1-4. 5+0MEGA++2.  >  >  + 

✓  1 . 5+(ALPHA2-.  5+OME6A++2.  >  +BIGD++2.  -BIGD 
6TWO=CO+ (BIGD+ (ALPHA1 +1 . 5+DMEGA++2. > 5> 

GTHREE=C6+  (BIGD+ (ALPHA 1+1 . 5+0MEGA++2. > -. 5> -2. ♦ (BIGD++2.  > 
♦DMEGA 

QUETW0=4.  +0MEGA++2.  -0MEGA0+*2.  -  (XI++2.  >  ♦  (GONE-GTWO> 

♦♦♦♦  COMPUTING  THE  COEFFICIENTS  OF  THE  CUBIC  FUNCTION  "C++ 2.“ 

COEFFC < 1 >  =GZER0++2 . 

COEFFC <2> =2. ♦GZERD+GUETWO 

COEFFC (3> =GUETW0++2.  +  <4. +B+OMEGA>  ++2. 

COEFFC.  <4>  =- 1 .  ♦  (XI ++4 .  )  +BETA2++2. 

•♦♦♦♦  CALLING  THE  SUBROUTINE  "ZPOLR"  TO  SOLVE  THE  CUBIC.  ♦♦♦♦♦♦ 

CALL  ZPOLR  (COEFFC » NDEG , ROOTS » I ER> 

■♦♦♦♦  ERROR  CHECK  ♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦ 

IF  (IER.NE. 0>  IERROR=IER 

♦♦♦♦♦  LOOP  TO  FIND  THE  REAL  ROOTS  OF  C++2  AND  COMPUTE  THE  MAX. 

♦♦♦♦♦  VALUE  OF  XOFTAU  FOR  EACH  OMEGA. 

DO  990  KK'=  1 »  3 

XMAX  (KK>  =0. 

CC  (KK>  =0. 

990  CONTINUE 
DO  99  J=1 *  3 

ROOTS I < J> =AI MAG (ROOTS ( J> > 

ROOTSR (  J>  =PEAL (ROOTS ( J> > 

IF  (ROOTS I <J> . EQ. 0.>  GO  TO  98 
GO  TO  99 
98  KSTOP=J+f 00 
KS TAP T*K STOP 

XI=B+SGRT (ROOTSR (J> ) ♦ (GTHREE+ (QUETWO+GZERO+ROOTSR (J> > 

✓  -4. ♦0MEGA+BETA2) ✓ (XI+BETA2) ++2. 

X2*-SG!RT  (ROOTSR  (J>  >  +BETR2+  (QUETWO+GZERD+RODTSR  (J>  >  ✓  (XI+BE 
TA2) ++2. 

PHI=ATAN2 (XI  a X2> 

PHASE=COS(PHI> 

TAU=0. 
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R=. 04/OMEGA 

DELTA=ATAN  <2.  ♦OMEGA^B-'  (0MEGA~2.  -OMEGA0^2.  >  > 

♦♦♦♦♦  LOOP  TO  CALCULATE  XOFTflU  ♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦ 


I 


q 


DO  88  K=KSTflRT*KSTOP 
THETfl=2. ♦OMEGA^TAU+PHI 

XOFTAU  (K)  =SQRT  (RDOTSR  <  J>  )  ♦CDS  (THETA)  +XI^  (  C0MEGA^2.  -OME& 
:  AQ^2.  >' 

✓  "  ♦♦2.+ (2. ♦B^OMEGA) ♦♦£.)♦♦-. 5^C0S(0MEGfl^ TAU+DELTfl>+ 

✓  Cl ♦XI#^2.  ♦ (C5+C£^CDS (2. ♦THETA) > ♦RODTSR (J) ♦ (C4^C0S  (OMEG 

A^ 

x  TAU+THETA)  +CO^COS  (OMEGA^TAU-THETA)  > ♦XI ♦SORT (RODTSR (J>  > 
TAU=TAU+R 
88  CONTINUE 

♦♦♦♦♦  LDDP  TO  FIND  THE  MAX.  OF  XDFTAU  ♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦ 

XMRX (J)  =XOFTAU  (KSTART) 

DO  7 7  L=KSTART » KSTOP 

I F  (XMRX  (  J)  .  LT.  XOFTflU  (L>  )  XMRX  (  J)  =XOFTAU  (L> 

77  CONTINUE 

XMRX  (J>  =XMAX  (J) ♦COS  (DELTA) 

CC  (  J)  =SQRT  (RDOTSR  (  J)  > 

CCC=CC(J) 

BJERK=AAZ^COS  (DELTA)  +CCC^C0^X2 

BJERK=BJERK>ARZ^ (XE^COS (DELTA) +X1 ♦SIN (DELTA) ) ♦CCC 
BJERK-X1 *BJERK 
XMRX ( J> *BJERK 
99  CONTINUE 
RETURN 
END 


THIS  IS  A  SUBROUTINE  TO  BE  USED  WITH  THE  PROGRAM  THAT 
COMPARES  OUR  RESULTS  WITH  THOSE  OF  PPOSREPETTI.  IT 
FINDS  THE  MAXIMUM  OSCILLATION  NEAR  RESONANCE. 


SUBROUTINE  RESNCE  (ALPHfll »  ALPHAE*  B>  OMEGAO*  XI  « W»  XMAX>  OMEGA) 
COMMON  AAZ 
COMMON/WR/AMPf COSPH 
REAL  XOFTflU (8 00) 

FAMP (X)  =  (OMEGA^B*. 1£5^SIN (2. ♦X) ♦ (XI/OMEGAO) ♦♦2.  )  ✓ <X1*D3+S1 
N(X))+  • 

s  SORT ( ( (OMEGA^B+. 1£5^SIN (£. ♦X) ♦ (XI^OMEGAO) ♦♦£. ) ✓ <XI*D3*SI 
• :  n  (x) ) ) 

/  ♦♦E.+1./D3) 

FPH1  (Y » X)  =Y^BO X-Dl  ^♦♦3.  +XI  ♦  (1 .  -D2+Y*Y>  ♦CDS  (X)  ♦.  25^Y^C0S  (2 

.  ♦X)  ♦ 

ei 
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/  CXI/CJMEGAO)  ♦♦2. 

UI=OMEGA/OME6AO 

♦♦♦♦  PRELIMINARY  EQUATIONS  ♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦ 

DEN=OMEGAO^-*2.  -4.  ♦OMEGA*^. 

D1  =  .  75^ ALPHA2- ALPHA 1  ♦  COME6AO*>*-2.  ♦  CALPHA1 -1 . 5^0MEGA^*2.  )  + 
.5^ 

/  CALPHA1 +1 . 5^0MEGA-^2. )  /DEN)  -1 . 5»0MEGA**2.  ♦  C.  25-  CALPHA1+1 
.  5* 

/  OMEGA^E.  )/DEN) 

D2=.  75^  CALPHA1-1 .  5*0MEGA»*2. )  /DEN+1 .  5^0MEGAO^-2.  ♦CALPHAl 

5^ 

✓  QMEGA^E. > 75 

D3=.  £5*  CALPHA1-7.  5*0MEGA*-*2. )  /HEN-*-.  5/OMEGA0^2.  ♦  CALPHAl -1 
.5^ 

✓  0MEGA^2.)-.25 

B0X=0MEGA^2.  -OME6AO»*2. -.  5*  CXI  /OMEGA 0>  ♦♦2.  ♦  C2.  ♦OMEGA**2. 
/  OMEGA0^2.)/DEN 

♦♦♦♦♦  PROGRAMMING  FOR  FINDING  THE  PHASE  ANGLE  AND  AMPLITUDE  ♦♦♦ 
♦♦♦♦♦  OF  OSCILLATION  ♦♦♦ 

J=0 

AAA=6. 282 

ACHNGE®. 06282 

FI  =FPHI  CFAMP  CAAA>  ..  AAA> 

1 0  AAA=AAA-ACHNGE 

IF  CAAA. LT. 0. >  GO  TO  40 
F2=FPHI  CFAMP < AAA) > AAA) 

SINCHK=F1^F2 

IF  CSINCHK.LT. 0.)  GO  TO  20 
15  J=0 

F1=F2 
GO  TO  10 

20  BBB=AAA+ACHNGE 
30  CCC=. 5^  CAAA+BBB) 

J=J+1 

IF  CJ.GT.500)  GO  TO  15 
.  FCHK=FPHI CFAMP  CCCC) > CCC) 

GUIT=. 0000001 

IF  CABS  CFCHK) . LE. QUIT)  GO  TO  40 
S I NCHK=F 1 ♦FCHK 
IF  CSINCHK.LT. 0.)  AAA=CCC 
IF  CSINCHK. GT. 0. )  BBB=CCC 
GO  TO  30 
40  PHI =CCC 


RMP=FRMP  <CCC> 

TRU=0. 

DELTfl=COS  CATANC  <2.  ♦B^OMEGA?  ✓  COMEGR^E.  -DME6R0^2.  >  >  > 

R=. 04/OMEGR 

♦♦♦♦  LOOP  TO  FIND  XOFTAU  ♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦ 

DD  200  M=1 f 2 

XDFTRU CM?  =RMP^CDS  COMEGA^TRU+PHI? ♦. 5^  CCDS (2. ♦DMEGRoTRU- 

✓  PHI>  ✓  C-l .  ♦DEN?  -COS  CPHI?  /DMEGR0~2.  ?  ♦XI^AMP-*-.  5^  C  C 

✓  RLPHR1  -1 . 5^DMEGA^2.  ?  'OMEGflO++2.  -  CRLPHfll  +1  -  5^0MEGA~2. 

>  ^  <— 1 .  ♦  i 

✓  DEN? ♦CDS  <2. ♦DMEGR^TRU+2. ♦PHI?  ? ♦AMP^2. 

TRU=TAU+R 

200  CONTINUE 

XMRX=XOFTRU  < 1 ? 

♦♦♦♦♦  LOOP  TD  FIND  THE  MAXIMUM  VALUE  DF  XDFTRU  ♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦♦ 
DD  300  M~1 »  2 

IF  CXMRX.LT. XDFTRU <N>?  XMRX=XDFTRU CN? 

300  CONTINUE 

XMRX=XMAX^DELTR 
COSPH=CDS  CPHI? 

XMRX=RMP^CDSPH 

RETURN 

END 
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